Density scaling as a property of strongly correlating viscous liquids 



OO 

o 
o 

(N 
S3 

S3 

(N 



o 



I 

S3 
O 

o 



> 
On 
OS 

(N 

rn 
O 
00 

O 



X 



Thomas B. Schr0der, Ulf R. Pedersen, and Jeppe C. Dyre 
DNRF Centre "Glass and Time, " IMFUFA, Department of Sciences, 
Roskilde University, Postbox 260, DK-4000 Roskilde, Denmark 
(Dated: June 2, 2008) 

We address a recent conjecture according to which the relaxation time r of a viscous liquid obeys 
density scaling (r = F(p 1 /T) where p is density) if the liquid is "strongly correlating," i.e., has 
almost 100% correlation between equilibrium virial and potential-energy fluctuations [Pedersen et 
al, PRL 100, 011201 (2008)]. Computer simulations of two model liquids - an asymmetric dumbbell 
model and the Lewis- Wahnstrom OTP model - confirm the conjecture and demonstrate that the 
scaling exponent 7 can be accurately predicted from equilibrium fluctuations. 

PACS numbers: 64.70.P- 



Undcrstanding the glass transition depends on under- 
standing the preceding highly viscous liquid phase, where 
the relaxation time upon cooling approaches and exceeds 
seconds [l[. Increasing the pressure also leads to much 
slower relaxations. The study of glass-forming liquids un- 
der high pressure has recently become popular, and many 
data are now available on the properties of the alpha and 
beta processes, etc., under pressure. If the density is p 
and T is temperature, the last few years have shown that 
for many highly viscous liquids the (alpha) relaxation 
time follows the scaling expression 



f(pVt) 



(1) 



The state of this rapidly developing field as of 2005 was 
summarized in the review Ref. [2| by Roland et al. that 
presented data for more than 50 liquids and polymers. 

Equation {T]) defines what is referred to as thermo- 
dynamic or density scaling. Density scaling is a recent 
discovery that, following pioneering works by Tolle and 
Dreyfus et al. [!, 0] , was proposed as a general princi- 
ple in 2004 in papers by Alba-Simionesco and co-workers 
and by Casalini and Roland [B|, Q • The former authors 
demonstrated data collapse at varying temperature and 
density following a more general expression than Eq. ([1]) 
used by the latter authors. 

Dreyfus and co-workers found 7 = 4 for ortho- 
tcrphcnyl (OTP) and argued that this could be linked 
to the r~ 12 repulsive term of the Lennard-Jones poten- 
tial. It turned out, however, that 7 = 4 is not a special 
exponent 0, leaving the question of the microscopic in- 
terpretation of 7 open. Coslovich and Roland very re- 
cently addressed this by computer simulations of binary 
Lennard-Jones like liquids where the exponent of the re- 
pulsive part of the potential took the values 8, 12, 24, 36 
[7( . These model systems obey density scaling and the 
exponent 7 is to a good approximation 1 /3 of the expo- 
nent characterizing the effective power law of the repul- 
sive core of the potential, as expected for an exact inverse 
power-law potential Q (see also Ref. [§]). 

Recently, simulations of the thermal equilibrium fluc- 
tuations of pressure, energy, and volume in different en- 



sembles revealed that these quantities correlate strongly 
for a number of model systems 0, [l(|. For instance, 
in the NVT ensemble (i.e., constant volume and tem- 
perature) the following systems are "strongly correlat- 
ing" in the sense that they show more than 90% corre- 
lation between virial (the non-kinetic part of the pres- 
sure) and potential energy: The standard Lennard-Jones 
liquid, a liquid with exponential short-range repulsion, 
the Kob- Andersen binary Lennard-Jones liquid, a seven- 
site united-atom model of toluene, and the model system 
studied below consisting of asymmetric "dumb-bell" type 
molecules. The correlations derive from the repulsive 
core of the intermolecular potential that, interestingly, 
dominate fluctuations even at zero and slightly negative 
pressure 1^, 11 1 . For the standard Lennard-Jones liquid 
the repulsive core is approximately described by a repul- 
sive r -18 term 1(3, 12 1. The exponent of the approximate 
power law depends weakly on state point. 

In view of the fact that both density scaling and the 
strong correlations reflects an effective inverse power-law 
of the repulsive core of the potential, the following con- 
jecture was proposed in Ref. [1(3]: A viscous liquid is 
strongly correlating if and only if it obeys density scal- 
ing to a good approximation. There is evidence indi- 
rectly supporting this: Highly viscous liquids with strong 
hydrogen bonds are not strongly correlating because of 
the "competing interactions" [lOj, and these liquids also 
do not obey density scaling very well @, [HI, 14 1. In 
the present publication simulations are presented that 
supports the conjecture and strengthens it by adding: 
For strongly correlating viscous liquids density scaling is 
obeyed with a scaling exponent 7 that can be determined 
from thermal equilibrium virial and potential energy fluc- 
tuations. This amounts to the simplest possible assump- 
tion; that it is the same part of the potential that domi- 
nates equilibrium fluctuations and flow dynamics. 

We performed NVT molecular dynamics simulations 
[15] ] of 512 asymmetric dumbbell molecules consisting of 
pairs of Lennard-Jones (LJ) spheres connected by rigid 
bonds. The dumbbells were parameterized to mimic 
toluene [3]. Charges of ±q (specified below) were ap- 
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q=0; p=1.166 g/cm A 3 
q=0; P = 2.0GPa 



O-Q q=0.5e; p= 1.006 g/cm A 3 




T[K] 

FIG. 1: Results from equilibrium molecular dynamics sim- 
ulations of 512 asymmetric dumbbell molecules with, respec- 
tively, a strong dipole moment (q = 0.5e, open symbols, 
three isochores), and zero dipole moment (q — Oe, filled sym- 
bols, three isochores and an isobar), (a) Correlation coeffi- 
cients, 7? = (AWAU) /y/{AW 2 ) {AU 2 ). (b) The 'slopes', 
7= ^{{AWY)/({AU) 2 ). 



plied to the LJ spheres. The model was simulated with 
two charges: i) q = corresponding to the simulations 
done in Refs. and [l(| • This version of the model is a 
"strongly correlating viscous liquid", ii) q — 0.5e (e be- 
ing the elementary charge) resulting in a dipole moment 
of 7.0D, i.e., almost 20 times stronger than in toluene, 
and almost 4 times stronger than water. The purpose of 
using such a large value of q was to break the correla- 
tions and thus to have a version of the model that is not 
strongly correlating. 

If the virial is denoted by W, and U is the potential 
energy, AW(t) = W(t) - (W) and AU(t) = U(t) - (U), 
where (...) indicates thermal average. The correlation co- 
efficient is denned by R = (AWAU) /^(AW 2 ) (AU 2 ). 
R is plotted for a several state points in Fig. [TJa). For 
q = (filled symbols) all investigated state points have 
R > 0.95. For q = 0.5e (open symbols) the correlation 
coefficient is significantly smaller; the Coulomb interac- 
tions do indeed break the correlations as expected ■ 

We define 7 = y/ ((AW) 2 ) / ((AU) 2 ). If R « 1 it fol- 
lows that AW(t) w jAU(t) in their instantaneous fluctu- 
ations [1(3], and consequently we refer to 7 as the 'slope'. 
According to our conjecture, the slope 7 is also the scaling 
exponent in Eq. (1). In Fig. [ljb) we show the slopes for 
all investigated state points. For q = (filled symbols) 
there is a small, but significant dependence on density 
and temperature. Thus if the conjecture is correct, den- 
sity scaling can not be exact: To apply density scaling a 
single value of 7 is needed, but the 7 we get from the flue- 
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FIG. 2: Arrhenius plot of the diffusion coefficient, D, for the 
asymmetric dumbbell model. 

tuations depends slightly on the state point In the 
following we consider for q = two values of the scaling 
exponent: the slope averaged over all state points with 
q = 0; 7 = 6.1, and the slope averaged over the three 
data sets with the smallest slopes for q = (p = 1.109 
g/cm 3 , p = 1.166 g/cm 3 and P=2.0GPa); 7 = 5.9, i.e., 
the 'best' compromise if we chose to ignore the p = 1.006 
g/cm 3 isochore. For q = 0.5e the slopes are less density 
dependent with a mean value 7 = 4.9. 

In the following we apply density scaling to the dif- 
fusion coefficient estimated from the long-time behavior 
of the mean-square displacem ent, (Ar 2 (i)), of the large 
spheres (the "phenyl group") [22j. The diffusion coeffi- 
cients for all state points studied are given in Fig. O 

Following Coslovich and Roland pj we apply den- 
sity scaling to the reduced diffusion coefficient, D* = 
(N/V) 1/3 (k B T/m)- 1 / 2 D where m is the mass of the 
molecules. In Fig. [3Ja) D* is plotted for q = as a func- 
tion of 1000p 7 /T with four different values of 7. Clearly, 
density scaling works neither with 7 = 7.0 or 7 = 5.0. 
Comparing the scaling with 7 = 6.1 to the data without 
scaling (filled symbols in Fig. [2|), we find here good data 
collapse; by far most of the density dependence is cap- 
tured by the density scaling with 7 = 6.1. With 7 = 5.9 
the data collapse is even better for three of the data sets, 
whereas one data set deviates slightly from the master 
curve comprised of these three sets. This is the isochore 
p = 1.006 g/cm 3 , i.e., the one that was ignored when 
choosing 7 = 5.9 (Fig. [2(b)). 

The conclusion drawn from Fig. [3] is two- fold: i) Den- 
sity scaling is approximate (as discussed above) - for a 
larger region of state points scaling will be less perfect, 
ii) For a given range of state points, the scaling exponent 
can be found by studying equilibrium fluctuations. 

In Fig. d] the reduced diffusion coefficients D* for 
q = 0.5e are plotted as a function of 1000p 7 /T with three 
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FIG. 3: The reduced diffusion coefficient, D* = 
{N/V) 1/3 (k B T/m)' 1/2 D, for the asymmetric dumbbell 
model with q = scaled according to Eq. (1), with four differ- 
ent scaling exponents: 7 = 7.0 (Upper set of curves, D* mul- 
tiplied by 100), 7 = 6.1 (Second set of curves, D* multiplied 
by 10), 7 = 5.9 (Third set of curves), and 7 = 5.0 (Lower set 
of curves, D* divided by 10). As a guide to the eye, the equa- 
tion D* = 4.07 x 10" 2 exp (-462/(T/p 5 ' 9 - 60.8)) is plotted 
as a fit to the three collapsing curves for 7 = 5.9. 



different values of 7. The value 7 = 4.9 chosen from the 
equilibrium fluctuations (Fig.QJb)), is found to be a rea- 
sonable scaling exponent. However, as conjectured, the 
data collapse achieved is inferior to that for the strongly 
correlating version of the model (q = 0) . 

To test the generality of our findings, we repeated 
the analysis for the Lewis- Wahnstrom OTP model (LW- 
OTP) [lElii (see also Refs. [H] and [HI). The results 
achieved for this model (Fig. [5]) are qualitatively similar 
to the results for the asymmetric dumbbell model with 
q = 0; (i) LW-OTP is strongly correlating (0.91 < R < 
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FIG. 4: The reduced diffusion coefficient, D* = 
(N/V) 1/3 (kBT/m)- 1/2 D, for the asymmetric dumbbell 
model with q = 0.5e scaled according to Eq. (1), with three 
different scaling exponents: 7 = 6.0 (Upper set of curves, D* 
multiplied by 10), 7 = 4.9 (Middle set of curves), and 7 = 4.0 
(Lower set of curves, D* divided by 10). 



0.92) [26(, (ii) the slope is slightly state point depen- 
dent, (iii) choosing the average slope as scaling exponent 
7 gives good data collapse. 

In summary, we have presented numerical evidence 
for the conjecture that density scaling is a property of 
strongly correlating viscous liquids. For the two strongly 
correlating models investigated density scaling applies 
with a scaling exponent that can be accurately predicted 
from the equilibrium fluctuations. This represents a step 
forward in the theoretical understanding of density scal- 
ing. In particular, the scaling exponent 7 should no 
longer be regarded as an empirical fitting parameter. In 
computer simulations 7 can be estimated directly from 
the equilibrium fluctuations. Via the fluctuation dissipa- 
tion theorem the equilibrium fluctuations determine the 
frequency-dependent linear thermoviscoelastic response 
functions - this fact provides a possible future route for 
independent experimental estimation of the scaling ex- 
ponent 7 dill- 

This work was supported by the Danish National Re- 
search Foundation's (DNRF) centre for viscous liquid dy- 
namics "Glass and Time." 
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FIG. 5: The reduced diffusion coefficient, D* = 
(N/V) 1/3 (k B T/rn)- 1/2 D 1 for LW-OTP @. For the zero 
pressure isobar the densities covers the interval 1.008 g/cm 3 
to 1.089 g/cm 3 Filled symbols: D* scaled according to Eq. (1) 
with 7 = 7.9 chosen from the slopes (inset). Open symbols: 
D* plotted without scaling (7 = 0). 
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